#######################################################

library(ggplot2)
library(dplyr)
library(optparse)

##########################################################################################
option_list <- list(
    make_option(c("--gobp_file"), type = "character"),
    make_option(c("--scriptPath"), type = "character"),
    make_option(c("--pct"), type = "character"),
    make_option(c("--type"), type = "character"),
    make_option(c("--logfc"), type = "character"),
    make_option(c("--out_path"), type = "character") 
)

if(1!=1){
    
    ## 通路富集的文件
    gobp_file <- "~/20231121_singleMuti/results/celltype_plot/diff_expression/one_vs_other.germ.pct_0.25.logfc_1.GOBP.tsv"

    ## 既往研究整理的代码
    scriptPath <- "~/20231121_singleMuti/scripts/scScalpChromatin"

    ## pct
    pct <- 0.25

    ## fc 
    logfc <- 1

    ## 输出
    out_path <- "~/20231121_singleMuti/results/celltype_plot/diff_expression/plot"
}

###########################################################################################
parseobj <- OptionParser(option_list=option_list, usage = "usage: Rscript %prog [options]")
opt <- parse_args(parseobj)
print(opt)

gobp_file <- opt$gobp_file
scriptPath <- opt$scriptPath
pct <- as.numeric(opt$pct)
logfc <- as.numeric(opt$logfc)
out_path <- opt$out_path
type <- opt$type

dir.create( out_path , recursive = T )

##########################################################################################
## 读入GOBP的文件
goRes_all <- read.delim( gobp_file , sep = "\t" )

##########################################################################################
## 已发表文献写好的脚本
source(paste0(scriptPath, "/plotting_config.R"))

##########################################################################################
## 画图的function
goPlot <- function(goRes = goRes , cell_type = cell_type , nterms = 10){
  cmap <- NULL
  border_color <- "black"
  barwidth <- 0.85
  enrichLimits <- c(0.0, 5.5)
  barLimits <- NULL

  # Plot GO results in bar plot form
  goRes$log2FoldEnrichment <- log2(goRes$Significant / goRes$Expected)
  goRes$log2FoldEnrichment <- ifelse(goRes$log2FoldEnrichment > enrichLimits[2], enrichLimits[2], goRes$log2FoldEnrichment)
  goRes$log2FoldEnrichment <- ifelse(goRes$log2FoldEnrichment < enrichLimits[1], enrichLimits[1], goRes$log2FoldEnrichment)
  goRes$threshPval <- ifelse(goRes$pvalue == "< 1e-30", 1e-30, as.numeric(goRes$pvalue))
  goRes$log10pval <- -log10(goRes$threshPval)
  if(is.null(cmap)){
    cmap <- cmaps_BOR$comet
  }
  if(!is.null(barLimits)){
    goRes$log10pval <- ifelse(goRes$log10pval < barLimits[2], goRes$log10pval, barLimits[2])
  }

  # Only plot the top nterms (reverse order to plot most significant at top)
  goRes <- goRes[1:nterms,]
  #goRes <- goRes[nrow(goRes):1,]
  #if(is.null(cmap)){cmap <- cmaps_BOR$comet}
  p <- ggplot(goRes, aes(y=log10pval, x=Term,fill=log2FoldEnrichment)) + 
    geom_bar(stat="identity", width=0.85,color="black")+xlim(rev(goRes$Term))+ 
    theme_test()+
    xlab("Term")+
    ylab("-log10pval")+
    ggtitle(cell_type)+
    theme (plot.title = element_text (size=12,hjust=0.5))+ 
    theme(axis.title = element_text(size = 14)) + 
    theme(axis.text = element_text(size = 11))+ 
    theme(plot.title = element_text(size = 17))+ 
    scale_fill_gradientn(colors=cmap, limits=enrichLimits) + 
    coord_flip()

  return(p)

}

## 每个细胞类型画图
for( cell_type in unique(goRes_all$cell_compare) ){
  
  ## 每个细胞类型分别画
  goRes <- subset( goRes_all , cell_compare == cell_type )

  ## 展示最显著的前10个通路
  p <- goPlot(goRes = goRes , cell_type = cell_type , nterms = 10)

  ## 打印图
  out_file <- paste0( out_path , "/" , "one_vs_other." , type , "." , "pct_" , pct, ".logfc_" , logfc , "." , cell_type , ".GOBP.pdf" )
  pdf(out_file , width = 10 , height = 5)
  print(p)
  dev.off()
}
